A "Tetris"-like model for the Compaction of Dry Granular Media 
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We propose a two-dimensional geometrical model, based on the concept of geometrical frustration, 
conceived for the study of compaction in granular media. The dynamics exhibits an interesting 
inverse logarithmic law that is well known from real experiments. Moreover we present a simple 
dynamical model of N planes exchanging particles with excluded volume problems, which allows to 
clarify the origin of the logarithmic relaxations and the stationary density distribution. A simple 
mapping allows us to cast this Tetris-like model in the form of an Ising-like spin systems with 
vacancies. 
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A granular system may be in a number of different mi- 
croscopic states at fixed macroscopic densities, and many 
unusual properties are linked to its non trivial packing 
UU ■ As pointed out by Edwards the role that the 
concept of free energy plays in standard thermal systems 
as Ising models, in granular media seems to be played by 
the "effective volume" , derived by a complex function of 
grain positions and orientations. In this way statistical 
mechanics provides theoretical concepts in the context of 
non-thermal systems. 

A recent experiment on the problem of density com- 
paction in a dry granular system under tapping, has 
shown H that density compaction follows an inverse log- 
arithmic law with the tapping number. 

Several approaches have been proposed to explain this 
behavior p|-|lO 12|, as geometrical models of "parking" 
pfj| , pl| or simple free- volume theories Q or the study of 
the dynamics of a frustrated lattice gas with quenched 
disorder subject to gravity and vibrations ||. 

In many seemingly different cases the logarithmic re- 
laxation proposed in B to describe experimental data is 
reproduced. Moreover the logarithmic law has turned out 
to be robust with respect to changes in the tapping pro- 
cedure H . This suggests that such a relaxation behavior 
is extremely general and not linked to specific properties 
of definite realizations. 

Here we introduce a purely geometrical model of sim- 
ple particles with several shapes on a lattice. We show 
that when subject to gravity and vibrations a logarithmic 
density relaxation [^| is found, due to the high entropic 
barriers (originated for geometrical reasons) to be passed 
by particles to improve global packing. 

We imagine a model, similar to the computer game 
Tetris in which neighboring grains can find different pack- 
ing volumes according to their relative geometrical ori- 
entations. Although one could imagine a rich variety of 
shapes and dimensions, like in the real computer game, 
it is useful, without loss of generality for the main fea- 



tures, to think just of a system of elongated particles 
which occupy the sites of square lattice tilted by 45° (see 
Fig. [I]), with periodic boundary conditions in the horizon- 
tal direction (cylindrical geometry) and a rigid plane at 
its bottom. In general the only interactions between the 
particles are the geometrical ones. Particles cannot over- 
lap and this condition produces very strong constraints 
(frustration) on their relative positions. For instance in 
the simplest case of two kind of elongated particles point- 
ing in two (orthogonal) directions, the frustration implies 
that two identical particles (pointing in the same direc- 
tion) cannot occupy neighboring sites in this direction. 
The particles are in principle allowed to rotate if at least 
three of their nearest neighbors are empty. This condi- 
tion is such that for sufficiently high densities the rotation 
events become negligible and the particles keep definitely 
their orientation. It is then reasonable, in the limit of a 
sufficiently large system, to consider an equal number of 
the two kind of particles. There is no other form of in- 
teraction between particles, and in this sense the model 
is purely geometrical. 

As stated, the particles are confined to a box and sub- 
jected to gravity. The effect of vibrations is introduced 
by allowing the possibility of moving also upwards, as 
explained below. 

The system is initialized by filling the container. The 
procedure of filling consists in inserting the grains at the 
top of the system, one at the time, and let them fall 
down, performing, under the effect of gravity, an oriented 
random walk on the lattice, until they reach a stable 
position, say a position in which they cannot fall further. 
This filling procedure is realized by the addition of one 
particle at the time and stops when no particles can enter 
the box from the top anymore. 

In our case the dynamics can be divided in two alter- 
nating steps. First, in a heating process (tapping) the sys- 
tem is perturbed by allowing the grains to move in any al- 
lowed directions with a probability p up to move upwards 
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(with < p U p < 0.5) and a probability pdown = 1 — Pup 
to move downwards. After each tapping has been com- 
pleted (i.e. a fixed number N moves per particle have 
been attempted with a fixed value of x = p U p/ Pdown) we 
allow the system to relax setting p up = 0. The relaxation 
process (Cooling) is supposed to be completed just when 
no particles can move anymore under just the effect of 
gravity, i.e. unless p up is switched on. After this relax- 
ation the system is in a stable static state and one starts 
again the cycle. We verified how the basic features of 
our model are very robust with respect to variations in 
the Monte-Carlo procedure. It is worth to stress how our 
dynamical procedure, is very close physically to the real 
processes of vibro-compaction [Q. Work is in progress 
to implement in our system the method proposed in fl^ j 
which allows for the simulation of a real tapping process. 

More precisely the single dynamical step consists of the 
following operations: 1) extracting with uniform proba- 
bility a grain; 2) extracting a possible movement for this 
grain among the 4 first neighbors (2 for the cooling pro- 
cess) according to the probabilities p up and Pdown] 3) 
move the grain if all the possible geometrical constraints 
with the neighbors are satisfied. 

We performed numerical simulations of the Tetris-like 
model in order to investigate its compaction properties. 
In particular we measured the density of the packing, i.e. 
the percentage of sites occupied with respect to the total 
number of sites, after each relaxation step and, in corre- 
spondence with real experiments, we plot the behavior of 
this density as a function of the number of taps. In order 
to avoid finite-size effects we considered systems with a 
linear size of at least L = 50 sites and, in order to be 
sure to observe bulk effects, we measured the density in 
the lower 25% of the system. 

Our main results on compaction are summarized in 
Fig. U which shows the evolution of the density, as a 
function of the number of taps, for different values of x 
and for a system of dimension L — 50. The different 
curves, obtained with a tap length of one iteration per 
particle, can be fitted according to the following inverse 
logarithmic law: 

Pitn) = P ~- 1 + B-to^Jr+l) (1) 

with = 1, a value Apoo = 0.25, which depends only 
on the loose packing density p to ~ 0.75, and two free 
parameters, B and a characteristic time t, for which we 
observe an algebraic dependence on x: 

t = Ax' 1 . (2) 

where 7 ~ 0.84 and A — 4.3. In this case r has the 
meaning of the minimum time over which one starts to 
observe a compaction process. Up to times t n << r, 
in fact, p(t n ) keeps practically the initial value. A com- 
plete and detailed analysis of these numerical results is 
reported in [flBI . 



Let us now briefly discuss how the system reaches the 
close-packing density which, just in the case of the sim- 
plest version with only two possible shapes, corresponds 
to a perfectly ordered state with unitary density. It is 
worth to stress again how this choice does not change the 
qualitative behavior of the system and an infinity of dis- 
ordered ground states can be obtained just allowing for a 
rich variety of shapes for the particles Jl5| . The approach 
to this state, realized by means of the two-step dynamics 
described above, represents a complex non-equilibrium 
process in which the system evolves alternatively with 
two different "temperatures" : a temperature T2 (heat- 
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ing process) such that e T 2 = 1 ^ ) p and a tempera- 
ture T\ = for the cooling process. The first step could 
be considered as a process going towards equilibrium in 
which detailed balance holds. Its features in many re- 
spects are very similar to the simple hard-square model 
||l6| . The step at zero temperature is an out of equilib- 
rium process which involves an irreversible positioning of 
the particles. Globally the microscopic reversibility and 
detailed balance are lost. 

In order to gain a deeper insight into the quoted loga- 
rithmic dynamical behaviors, let us introduce and discuss 
a simple model which describes the evolution of a system 
of particles which hop on a lattice of n = 0, N stacked 
planes according to the ideas of "parking" introduced in 
|10| , ^T| . We consider a system of particles which can move 
up or down between N layers in such a way that their 
total number is conserved. We ignore the correlations 
among particles rearrangements and the problem related 
to the mechanical stability of the system. The master 
equation for the density on a generic plane n, except for 
the n — plane, is given by: 

dtPn = (1 - Pn)D(p n )[p n -i ■ p up + p n+ i ■ Pdown)} + 

-Pn[(l - Pn-l)D(p n -l)Pdown + (3) 
+ (1 - Pn+l)D(p n+ i)p U p} 

where Pdown and p up have been defined above for the 
Tetris-lihe model. D(p n ) is a sort of mobility for the 
particles given by the probability that the particle could 
find enough space to move. Apart form other effects 
it takes mainly onto account the geometrical effects of 
frustration, i.e. the fact that the packing prevents the 
free move of the particles. In a naive way one could 
imagine a functional form like D(p n ) — p n (l — p n ') ob- 
tained by considering only the nearest neighbors interac- 
tions in the Teim-like model. It is easy to realize that 
such an approach does not account for the complexity 
of the problem where the packing at high densities cre- 
ates long range correlations in the system, and, using 
this functional form, the equations show a trivial expo- 
nential relaxation. Non trivial results are obtained with 
a careful choice of the functional form for D(p n ) which 
takes onto account the cooperative effects on the dynam- 
ics generated by the frustration. This functional form 
can be obtained by evaluating the rate specifying how 



2 



many steps are needed in a frustrated system, i.e. the 
Tetris-\ike model, with respect to a non-frustrated one, 
to achieve a rearrangement in a new configuration. Wc 
do not report here the complete calculation of this rate 
which involves the evaluation of the configurational en- 
tropy of the frustrated and of the non-frustrated systems 
respectively, and we refer to [Q. It is nevertheless par- 
ticularly enlightening to consider the example N = 2. In 
this case the system can be reduced to a one-dimensional 
chain in which the average length of the filled intervals 
turns out to be given by < I >~ p/(l — p). The number 
of steps to move a certain particle, related to the number 
of steps necessary to move the entire interval (of average 
length < I >) aside that particle, will be the order of 
N p ~ exp < I >. The general form of D(p n ), although 
very complicate, must then include a term like 

D(p n ) = D exp[-p n /(l - p n )]. (4) 

We checked this functional form, which can be seen how 
the outcome of a free- volume theory for granular media 
p0j , p~2|JTr| , by means of specific simulations for that quan- 
tity on the TetrisASke model Jig] , 

In the general case with arbitrary N we obtained the 
exact asymptotic stationary solution for the density on 
each plane. This solution is in an implicit form and for 
that we refer to jl5| . It is possible nevertheless to extract 
the approximate explicit behaviors. In particular, if M is 
the total "mass" of the system, i.e. the maximal number 
of planes which can be completely filled, one gets: 

( p^ ~ 1 - 1/[(M - fc) • log(l/x)] for k«M 

I (5) 

[ p™ ~ e (M-fc)-log(l/a:) for k M 

The stationary solution tends thus to a step function 
9(k - M) in the limit x -> 0. 

Let us now comment on the dynamical behavior of the 
system, i.e. the relaxation towards the stationary solu- 
tion. We start by considering the simplest case with just 
N = 2 planes. In the limit x = p U p/Pdown << 1 and 
for a sufficiently high total density, p = 1 — e (e << 1), 
one can easily prove |L5| that the asymptotic equilibrium 
lower plane density behaves like ~ 1 + l/log[x/(e)] 
where /(e) = 2e/(l - 2e) exp[-(l - 2e)/2e]. The dynam- 
ical equation for pi can be written exactly. In the limit 
px >> (1 — 2e) (which holds for sufficiently long times) 
this equation exhibits a very simple form as: 

dtpi = B(x, e)(l - pi)D( Pl ) - A(x, e)pi (6) 

mthB(x,e) = (l-2e)/(l+x) and A(x,e) = 2eexp[-(l- 
2e)/2e]x/ (1 + x). This equation has the same form of the 
one dimensional "parking problem" studied in [jl0| whose 
absorption and desorption parameters are now written in 
terms of the global density of the system, 1 — e, and of 
the vibration amplitude ratio x. It exhibits a logarith- 
mic solution up to times of the order of to ~ 1/A(x, e) 
Q. Later on, when the density approaches its steady 



state value, the first term on the right-hand side of eq.(|6|) 
becomes negligible with respect to the loss term and an 
exponential saturation becomes dominant. Here we just 
note that to grows inversely proportional to x but has 
an essential singularity for e going to zero. So for low 
enough amplitude vibrations or high enough densities 
the logarithmic region extends actually up to any ex- 
perimentally observable time. The cooperative effect of 
interaction among the different planes makes the times 
over which one observes the logarithmic relaxation longer 
and longer. Crucial for this effect is the value of e, i.e. 
the asymptotic difference of density between two adja- 
cent planes. In the general case of N planes, one has from 
eq.(g) that in the bulk, for finite values of x, tk — x M ~ k , 
i.e. the are exponentially small in M. We then expect 
that the logarithmic relaxations extends up to times of 
the order of x~ M (see Jl5| for a detailed discussion of this 
point). 

Let us now notice a further aspect of our model. The 
two-step dynamics of our model may be easily interpreted 
in terms of a Glauber dynamics for an Hamiltonian with 
Ising-like variables. In this language the geometrical 
model is mapped into the following Ising-like Hamilto- 
nian with vacancies in the limit J — > oo: 

H = }^J(SjSj - dij(Si + Sj) - VjUiUj +g^y(i) (7) 

where n, =0,1 are occupancy variables and Si = ±1 
are spin variables that corresponds to the twofold orien- 
tation of the particles, = ±1 are fixed non-random 
bond fields with an ordered structures: = 1 for bonds 
along one direction of the lattice and = — 1 for bonds 
in the other. In the gravitational term, ~9^2 i y{i) 1 g is 
the gravity, and y(i) is the ordinate of the lattice site 
i. It is easy to realize that the sum of the ay converg- 
ing on each single site is zero. This implies that the 
ground state of Hamiltonian (^) is perfectly antiferro- 
magnetic if the densities of the two kinds of particles are 
equal. This state is reached just when all the sites of the 
lattice are occupied, so = 1 Vi. This mapping, and 
the ones for the models with a variety of shapes which 
lead to Potts-like Hamiltonians, are particularly useful 
as starting points for an analysis of these systems in a 
thermodynamic framework |lq| . 

In this paper we have introduced a very simple geo- 
metrical model in order to describe the phenomenon of 
compaction in dry granular media. It takes in account 
excluded volume effects, say the geometrical constraints 
which are felt by granular media during the relaxation to- 
wards the highest density optimal packing configuration. 
When subjected to Monte Carlo vibrations, defined by 
a diffusive dynamics, it exhibits a density compaction 
after tapping which reproduces the inverse logarithmic 
behavior found in both experiments |5j and other models 

This Tetris-\ike model can be easily generalized by in- 
troducing an arbitrary fixed number of shapes for the 
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particles which correspond to complicated matrices for 
the particle-particle interactions. This kind of general- 
ization does not change the qualitative structure of the 
relaxation but it could account for other effects of disor- 
der in granular media: segregation, hysteresis etc. Eaj. 

Furthermore we presented a simple dynamical model 
of N planes exchanging particles with excluded volume 
effects. For this model we have found the exact stationary 
density distribution and we have shown how, without 
loss of generality with respect to the choice of particular 
geometrical constraints, it allows for an explanation of 
the inverse logarithmic law for compaction. 

It is moreover interesting that the pure geometrical 
model presented here can be mapped into a simple Hamil- 
tonian formalism of an Ising antiferromagnet. This con- 
nects our work to previous works [p|,p|Jl0|] introduced to 
discuss different aspects of granular media phenomenol- 
ogy, and could open the way to their systematic anal- 
ysis [fL5[ , Acknowledgments We are indebted with 
A. Coniglio for useful suggestions. We thank P.G. De 
Gennes for bringing to our attention his preprint. 
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FIG. 1. Schematic picture of one particular configuration of 
the grains considered in the simplest version of the Tefris-like 
model. The two types of particles have to fulfill only geomet- 
rical constraints in their dynamics. As shown in the figure 
these constraints are due to the impossibility for the particles 
to overlap. 




FIG. 2. Logarithmic behavior of the density of the pack- 
ing, measured in the lower 25% of the system, as a function of 
tapping number t n , for five different values of amplitude vi- 
brations x = pup/pdomn = 0.001, 0.01, 0.03, 0.1, 0.5, from bot- 
tom to top. The superimposed logarithmic fit curves, given 
by eq.(0), were proposed to describe experimental data. 
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